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Abstract. Massive stars influence the surrounding universe far out of proportion to their num- 
bers through ionizing radiation, supernova explosions, and heavy element production. Their 
formation requires the collapse of massive interstellar gas clouds with very high accretion rates. 
We discuss results from the first three-dimensional simulations of the gravitational collapse of 
a massive, rotating molecular cloud core that include heating by both non-ionizing and ioniz- 
ing radiation. Local gravitational instabilities in the accretion flow lead to the build-up of a 
small cluster of stars. These lower-mass companions subsequently compete with the high-mass 
star for the same common gas reservoir and limit its overall mass growth. This process is called 
fragmentation-induced starvation, and explains why massive stars are usually found as members 
of high-order stellar systems. These simulations also show that the H n regions forming around 
massive stars are initially trapped by the infalling gas, but soon begin to fluctuate rapidly. Over 
time, the same ultracompact H II region can expand anisotropically, contract again, and take 
on any of the observed morphological classes. The total lifetime of H II regions is given by the 
global accretion timescale, rather than their short internal sound-crossing time. This solves the 
so-called lifetime problem of ultracompact H II region. We conclude that the the most significant 
differences between the formation of low- mass and high-mass stars are all explained as the result 
of rapid accretion within a dense, gravitationally unstable flow. 
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1. Introduction 



High-mass stars form in denser and more massive cloud cores (Motte et al. 2008) 



than their low- mass counterparts (Myers et al.||1986). High densities result in the large 



accretion rates, exceeding 10 -4 M Q yr _1 , required for massive stars to reach their final 
mass before exhausting their nuclear fuel ( Keto fc Wood|2006 ). High densities also result 
in local gravitational instabilities in the accretion flow, resulting in the formation of 



multiple additional stars (Klessen & Burkert 2000, 2001 Klessen 2001 Kratter & Matzner 
2006). Young massive stars almost always have companions ( Ho fc Haschick||1981 ), and 
the number of their companions significantly exceeds those of low-mass stars ( Zinnecker 



& Yorke 2007). Such companions influence subsequent accretion onto the initial star 
(Krumholz et al. 2009). Observations show an upper mass limit of about 100 M . It 



remains unclear whether limits on internal stability or termination of accretion by stellar 
feedback determines the value of the upper mass limit ( Zinnecker fc Yorke||2007 1 . 
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H II regions form around accreting protostars once they exceed ~ 10 M Q , equivalent 
to a spectral type of early B. Thus, accretion and ionization must occur together in the 
formation of massive stars. The pressure of the 10 4 K ionized gas far exceeds that in the 
10 2 K accreting molecular gas, creating unique feedback effects such as ionized outflows 
( |KeTo][2OT2l[2503l[2007 ). 

Around the most luminous stars the outward radiation pressure force can equal the 
inward gravitational attraction. A spherically symmetric calculation of radiation pressure 
on dust yields equality at just under 10 M (Wolfire & Cassinelli 1987[ ). However, the 
dust opacity is wavelength dependent, the accretion is non-spherical, the mass-luminosity 
ratio is different for multiple companions than for a single star, and to stop accretion more 
than static force balance is required. The momentum of any part of the accretion flow 



must be reversed (Larson & Starrfield 1971 



et al. 1995). Observations by Keto & Wood 



Kahn 1974 Yorke & Krugel||1977 Nakano 



( 2006 ) provide evidence for the presence of 



all these mitigating factors, and numerical experiments combining some of these effects 



(Yorke & Sonnhalter 2002 Krumholz ct al. 2007) confirm their effectiveness, showing 



that radiation pressure is not dynamically significant below the Eddington limit. 

The most significant differences between massive star formation and low-mass star 
formation seem to be the clustered nature of star formation in dense accretion flows and 
the ionization of these flows. We present results from three-dimensional simulations by 
Peters et al. ( 2010a|b|c ) of the collapse of molecular cloud cores to form a cluster of 
massive stars that include ionization feedback. These calculations are the first ones that 
allow us to study these effects simultaneously. 



2. Modeling High-Mass Star Formation 

Our discussion is based on a series of recent numerical simulations by |Peters et al.| 



( 2010a|b|c ). They are based on a modified version of the adaptive-mesh code FLASH 
(Fryxell et al. 2000) that has been extended to include sink particles representing pro- 



tostars (Federrath et al. 20101. The protostars evolve following a prestellar model that 



determines their stellar and accretion luminosities as function of protostellar mass and 
accretion rate. We set the stellar luminosity with the zero-age main sequence model by 



Paxton (20041 and the accretion luminosity by using the tables by Hosokawa & Omukai 
( 2009 1 . The ionizing and non- ionizing radiation from the protostars is propagated through 



the gas using an improved version of the hybrid characteristics raytracing method on the 



adaptive mesh developed by Rijkhorst et al. (2006). In some calculations, secondary sink 



formation is suppressed with a density-dependent temperature floor to prevent runaway 



collapse of dense blobs of gas. For further details, consult Peters et al. (2010a). 

The simulations start with a 1000 M molecular cloud. The cloud has a constant 
density core of p — 1.27 x 10 -20 g cm~ 3 within a radius of r = 0.5 pc and then falls off as 
r -3/2 un tji r — 1.6 pc. The initial temperature of the cloud is T = 30 K. The whole cloud 
is set up in solid body rotation with an angular velocity tu = 1.5 x 10 -14 s _1 corresponding 
to a ratio of rotational to gravitational energy j3 = 0.05 and a mean specific angular 
momentum of j = 1.27 x 10 23 cm 2 s _1 . Peters et al. ( 2010a|b[c ) follow the gravitational 
collapse of the molecular cloud with the adaptive mesh until they reach a cell size of 
98 AU. Then sink particles are created at a cut-off density of p cr i t = 7x 10~ 16 gcm~ 3 . 
All gas within the accretion radius of r S j n k = 590 AU above p cr it is accreted to the sink 
particle if it is gravitationally bound to it. The Jeans mass on the highest refinement level 
is Mj H a,ns = 0.13 Mq. A summary of the three simulations discussed in this proceedings 



jeans 



article is provided in Table [TJ 
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Table 1. Overview of collapse simulations. 



Name Resolution Radiative Feedback Multiple Sinks M s i n k s (Mq) iV s i n ks M max (Mg) 



Run A 98 AU 
Run B 98 AU 
Run D 98 AU 



yes 
yes 
no 



no 

yes 

yes 



72.13 


1 


72.13 


125.56 


25 


23.39 


151.43 


37 


14.64 



3. Fragmentation-Induced Starvation 

In this section we compare the protostellar mass growth rates from three different 
calculations. Run A only allows for the formation of a single sink particle, Run B has 
multiple sinks and radiative heating, while Run D has multiple sinks but no radiative 



heating. As already discussed by Peters et al. (2010a), when only the central sink particle 
is allowed to form (Run A), nothing stops the accretion flow to the center. Figure fl] shows 
that the central protostar grows at a rate M « 5.9 x lO _4 M yr _1 until we stop the 
calculation when the star has reached 72 M . The growing star ionizes the surrounding 
gas, raising it to high pressure. However this hot bubble soon breaks out above and below 
the disk plane, without affecting the gas flow in the disk midplane much. In particular, it 
cannot halt the accretion onto the central star. Similar findings have also been reported 
from simulations focussing on the effects of non-ionizing radiation acting on smaller scales 
( |Yorke fc Sonnhalter|[2002| |Krumholz et al.|2007| |2009| |Sigalotti et"al~1|2009[ ) . Radiation 
pressure cannot stop accretion onto massive stars and is dynamically unimportant, except 
maybe in the centers of dense star clusters near the Galactic center. 

The situation is different when the disk can fragment and form multiple sink particles. 
Initially the mass growth of the central protostar in Runs B and D is comparable to the 
one in Run A. However, as soon as further protostars form in the gravitationally unstable 
disk, they begin to compete with the central object for accretion of disk material. Unlike 



in the classical competitive accretion picture (Bonnell et al. 2001a, 2004), it is not the 



most massive object that dominates and grows disproportionately fast. On the contrary, 
it is the successive formation of a number of low-mass objects in the disk at increasing 
radii that limits subsequent growth of the more massive objects in the inner disk. Material 
that moves inwards through the disk due to viscous and gravitational torques accretes 
preferentially onto the sinks at larger radii. 

Similar behavior is found in models of low-mass protobinary disks, where again the 
secondary accretes at a higher rate than the primary. Its orbit around the common center 
of gravity scans larger radii and hence it encounters material that moves inwards through 
the disk before the primary star. This drives the system towards equal masses and circular 



orbits (Bate 2000). In our simulations, after a certain transition period hardly any gas 



makes it all the way to the center and the accretion rate of the first sink particle drops 
to almost zero. This is the essence of the fragmentation-induced starvation process. In 
Run B, it prevents any star from reaching a mass larger than 25 Mq. The Jeans mass in 
Run D is smaller than in Run B because of the lack of accretion heating, and consequently 
the highest mass star in Run D grows to less than 15 M Q . 

Inspection of Figure [T] reveals additional aspects of the process. We see that the total 
mass of the sink particle system increases at a faster rate in the multiple sink simulations, 
Runs B and D, than in the single sink case, Run A. This is understandable, because as 
more and more gas falls onto the disk it becomes more and more unstable to fragmen- 
tation, so as time goes by additional sink particles form at larger and larger radii. Star 
formation occurs in a larger volume of the disk, and mass growth is not limited by the 
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Figure 1. left Total accretion history of all sink particles combined forming in Runs A, B, and 
D. right Instantaneous accretion rate as function of time of the first sink particle to form in the 
three runs. 

disk's ability to transport matter to its center by gravitational or viscous torques. As a 
result the overall star- format ion rate is larger than in Run A. 

Since the accretion heating raises the Jeans mass and length in Run B, the total number 
of sink particles is higher in Run D than in B, and the stars in Run D generally reach 
a lower mass than in Run B. These two effects cancel out to lead to the same overall 
star formation rate for some time. Eventually, however, the total accretion rate of Run B 
drops below that of Run D. At time t rj 0.68 Myr the accretion flow around the most 
massive star has attenuated below the value required to trap the H n region. It is able 
to break out and affect a significant fraction of the disk area. A comparison with the 
mass growth of Run D clearly shows that there is still enough gas available to continue 
constant cluster growth for another 50 kyr or longer, but the gas can no longer collapse 
in Run B. Instead, it is swept up in a shell surrounding the expanding H n region. The 
figure furthermore demonstrates that, although the accretion rates of the most massive 
stars (M > 10 M Q ) steadily decrease, the low-mass stars (M < 10 M ), which do not 
produce any significant H n regions, keep accreting at the same rate. 



4. Properties of the H n Regions 

In all calculations by Peters et al. (2010a|c), the H II regions are gravitationally trapped 



in the disk plane but drive a bipolar outflow perpendicular to the disk. The highly 
variable rate of accretion onto protostars as they pass through dense filaments causes fast 
ionization and recombination of large parts of the interior of the perpendicular outflow. 
The H II regions around the massive protostars do not uniformly expand, but instead 
rapidly fluctuate in size, shape and luminosity. 

We can directly compare these numerical models with radio observations of free-free 
continuum, hydrogen recombination lines, and NH3(3,3) rotational lines by generating 
synthetic maps ( Peters et al. 2010a|b I . The simulated observations of radio continuum 
emission reproduce the morphologies reported in surveys of ultracompact H II regions 
( Wood &: Churchwell|1989 Kurtz et al.|1994 ). Figure [2] shows typical images from Run B 
to illustrate this point. It is important to note that even the correct relative numbers of 
the different morphological types are obtained. Table [2] shows the morphology statistics 



of UC H II regions in the surveys of Wood & Churchwell (1989) and Kurtz et al. (1994) 



as well as from a random evolutionary sample from Runs A and B of 500 images for each 
simulation. While the statistics of the cluster simulation B agrees with the observational 
data, this is not the case for Run A, in which only one massive star forms. For further 



discussion, see Peters et al. (2010b) 
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Figure 2. H II region morphologies. This figure shows ultracompact H II regions around massive 
protostars in Run B at different time steps and from different viewpoints. The cluster is assumed 
to be 2.65 kpc away, the full width at half maximum of the beam is 0'.'14 and the noise level is 
10 _3 Jy. This corresponds to typical VLA parameters at a wavelength of 2 cm. The protostellar 
mass of the central star which powers the H II region is given in the images. The H II region 
morphology is high l y variable in time and shape, taking the form of any observed type ( Wood 
fc Churchwell|1989 Kurtz et al.|1994 1 during the cluster evolution. 



Table 2. Percentage Frequency Distribution of Morphologies 



Type 



Wood & Churchwell (19891 



Kurtz et al. 



(1994) Run A Run B 



Spherical/Unresolved 
Cometary 
Core-halo 
Shell-like 
Irregular 



43 
20 
16 

4 

17 



55 


19 


60 ± 5 


16 


7 


10 ± 5 


9 


15 


4 ± 2 


1 


3 


5 ± 1 


19 


57 


21 ± 5 



The H II regions in the model fluctuate rapidly between different shapes while accretion 
onto the protostar continues. When the gas reservoir around the two most massive stars 
is exhausted, their H n regions merge into a compact H II region, the type that generally 



accompanies observed ultracompact H II regions (Kim & Koo 2001 ). These results suggest 
that the lifetime problem of ultracompact H II regions ( Wood fc Churchwell|1989 ) is only 
apparent. Since H n regions embedded in accretion flows are continuously fed, and since 
they flicker with variations in the flow rate, their size does not depend on their age until 
late in their lifetimes. 

To compare more directly to observations of the time variability of H II regions, [Peters] 
et al. (2010b) analyze a few time intervals of interest at a resolution of 10 years. They 



find that when the accretion rate to the star powering the H II region has a large, sudden 
increase, the ionized region shrinks, and then slowly re-expands. This agrees with the 
contraction, changes in shape, or anisotropic expansion observed in radio continuum 
observations of ultracompact H II regions over intervals of ~ 10 yr (|Franco-Hernandez fc 
Rodriguez||2"0"04| [Rodriguez et al.||2007[ |Galvan-Madrid et al.||2008| ). Figure |3| shows the 
2-cm continuum flux, the characteristic size of the H II region, and the rate of accretion 
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Figure 3. Time variability of the accretion flow and H II region. The image shows accretion rate 
onto the sink particle M s j n k (sink particle, in units of 10 -6 Mq yr" 1 ), corresponding diameter 
H of the H II region (in units of 10 AU), as well as the resulting 2-cm continuum flux F2 cm (in 
units of 10 mjy) for two selected time periods in Run B. 



onto the star. In the left panel, the H II region is initially relatively large, and accretion 
is almost shut off. A sudden strong accretion event causes the H n region to shrink and 
decrease in flux. The star at this moment has a mass of 19.8 M Q . In the right panel, 
the star has a larger mass (23.3 M©), the H n region is initially smaller, and the star 
is constantly accreting gas. The ionizing-photon flux appears to be able to ionize the 
infalling gas stably, until a peak in the accretion rate by a factor of 3 and the subsequent 
continuos accretion of gas makes the H II region to shrink and decrease in flux. The 
H II region does not shrink immediately after the accretion peak because the increase 
is relatively mild and the geometry of the infalling gas permitted ionizing photons to 
escape in one direction. These results show that observations of large, fast changes in 
ultracompact H n regions ( Franco-Hernandez fe Rodriguez||2004 Rodriguez et al.||2007 



Galvan-Madrid et al.|2008) are controlled by the accretion process 



5. Comparison with W51e2 



The Peters et al. (2010a) model can be compared with the well-studied ultracompact 
H II region W51e2 ( |Zhang et al]|1998[ |Keto fc Klaassen||2008"l ). In Figure [I] we show 
simulated and observed maps of the NH3(3, 3) emission, the 1.3 cm thermal continuum 
emission, and the H53a radio recombination line. The simulated maps were made at a 
time when the first star in Run B has reached a mass of 20 M© . 

The brightest NH 3 (3, 3) emission reveals the dense accretion disk surrounding the most 
massive star in the model, one of several within the larger-scale rotationally flattened 
flow. The disk shows the signature of rotation, a gradient from redshifted to blueshiftcd 
velocities across the star. A rotating accretion flow is identifiable in the observations, 
oriented from the SE (red velocities) to the NW (blue velocities) at a projection angle of 
135° east of north (counterclockwise). 

The 1.3 cm radio continuum traces the ionized gas, which in the model expands per- 
pendicularly to the accretion disk down the steepest density gradient. As a result, the 
simulated map shows the brightest radio continuum emission just off the mid-plane of 
the accretion disk, offset from the central star rather than surrounding it spherically. 
Continuum emission in the observations is indeed offset from the accretion disk traced in 
ammonia. The NH 3 (3,3) in front of the H n region is seen in absorption and red-shifted 
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Figure 4. Comparison of line and continuum emission simulated from the model (upper panels) 
and actually observed from the W51e2 region (lower panels). The left panels show the NH3(3,3) 
line emission strength in white contours, the molecular line velocities as the background color, 
and the 1.3 cm free-free continuum from ionized gas in red contours. The right panels show the 
H53 q recombination lin e velocities from the ionized gas. The molecular lin e obs ervations are 



from [Zhang ct a lT| (]1998|) and t he H53a observations from Keto & Klaasscn (2008). For further 
details see|Peters et al.|(|2010a[) . 



by its inward flow toward the protostar. The density gradient in the ionized flow deter- 
mines the apparent size of the H II region. Therefore the accretion time scale determines 
the age of the H II region rather than the much shorter sound-crossing time. 

Photoevaporation of the actively accreting disk supplies the ionized flow. Therefore, the 
ionized gas rotates as it flows outward, tracing a spiral. An observation that only partially 
resolves the spatial structure of the ionized flow sees a velocity gradient oriented in a 
direction between that of rotation and of outflow, as shown in the simulated observation. 
The observed H53a recombination line ( Keto fc Klaassen|2008 1 in Figure 4 indeed shows 
a velocity gradient oriented between the directions of rotation and the outflow. 



6. Conclusions 

Numerical simulations of high-mass star formation regions are now able to resolve the 
collapse of massive molecular cloud cores and the accretion flow onto the central group 
of protostars while at the same time treating the transport of ionizing and non-ionizing 
radiation. This allows us to consistently follow the dynamical evolution of the H n regions 
that ubiquitously accompany the birth of massive stars. We find that the accretion flow 
becomes gravitationally unstable and fragments. Secondary star formation sets in and 
consumes material that would otherwise be accreted by the massive star in the center. We 
call this process fragmentation-induced starvation. It determines the upper mass limit 
of the stars in the system and explains why massive stars are usually found as members 
of larger clusters. These simulations furthermore show that the H II regions forming 
around massive stars are initially trapped by the infalling gas. But soon, they begin 
to fluctuate. Over time, the same ultracompact H n region can expand anisotropically, 
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contract again, and take on any of the observed morphological classes. The total lifetime 
of H II regions thus is given by the global accretion timescale, rather than their short 
internal sound-crossing time. This solves the so-called lifetime problem of ultracompact 
H II regions. 
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